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o : 

We investigate the dynamics of a three-state stochastic lattice gas, consisting of holes and two 
■ oppositely "charged" species of particles, under the influence of an "electric" field, at zero total 

charge. Interacting only through an excluded volume constraint, particles can hop to nearest neigh- 
bour empty sites. Using a combination of Langevin equations and Monte Carlo simulations, we study 
the probability distributions of the steady-state structure factors in the disordered phase where homo- 
geneous configurations are stable against small harmonic perturbations. Simulation and theoretical 
results are in excellent agreement, both showing universal asymmetric exponential distributions. 
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Analyzing structure factors and correlations in many-particle systems is a standard way to study collective behaviour 
in real experiments, computer simulations and theoretical frameworks. For example, for a system with short range 
microscopic interactions, placed in thermal equilibrium, there are no long range spatial correlations in general. Their 
presence is a typical signal that the system is at a critical point. On the other hand, when such systems are driven 
' into non-equilibrium steady states, long range correlations are often observed jjj . Another good example is the driven 
7— I , Ising lattice gas [^]j3|. The role of the drive ("external field") is to bias hopping rates along a particular direction 
on the lattice. In addition to many other unexpected features such as non-Hamiltonian fixed points controlling the 
critical behavior Q, the system exhibits long range spatial correlations at all temperatures above criticality, as a 
result of the breakdown of the traditional fluctuation-dissipation relations [^J. In momentum space, this appears as a 
discontinuity singularity of the structure factor at the origin || . 

In this Paper, we revisit a simple, non-equilibrium model |?| displaying similar features and examine the steady state 
structure factors from a more general point of view. In particular, we will construct the full distribution of various 
structure factors, using both Monte Carlo simulations and Langevin equations following the method in ref The 
model can be considered as a generalization of the KLS model, in a fashion similar to an equilibrium Ising model 
being extended to spin-1 Q or Potts |l^] models, in that it involves two species of particles (and holes). Referred to 
as +'s and — ! s, they are driven in opposite directions, subject to periodic boundary conditions. To keep the model 
simple, we suppress the usual Ising nearest neighbor interaction and retain only the excluded volume constraint and 
the bias. In the simplest scenario, charge exchange (+ <-» — ) is not allowed. Phase transitions in this multi-species 
model, in both one and two dimensions, have been investigated extensively [fTT| ^6|. Monte Carlo simulations [ll| in 
two dimensions and mean-field studies |T2|jr^ ] show that there is a transition, controlled by particle density and drive, 
from a spatially homogeneous (disordered) phase to a charge segregated one, where the excluded volume constraint 
leads to the mutual blocking of particles. We will focus on the disordered phase of the system, where we have a, sound 
analytic understanding of the dynamics in terms of Langevin equations. Performing Monte Carlo simulations and 
measuring the full distributions of several structure factors, we find that the Gaussian approximation is extremely 
effective in predicting these properties of our non-equilibrium system. In contrast to the earlier study ]§| , the structure 
factors here form a 2 x 2 matrix, due to the presence of two species of particles, and the real and imaginary parts of 
the off-diagonal elements exhibit a characteristic asymmetry. The measurement of these distributions can be used to 
extract considerable information about the structure of the underlying linearized equations of motion, independent of 
model details. 
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II. THE MICROSCOPIC MODEL AND THE SIMULATIONS 



We consider a two dimensional fully periodic lattice with L x L sites, each of which can be empty or occupied 
by a single particle. Since there are two species we need two occupation variables n+ and n~, with n being or 
1, depending on whether a positive or negative particle is present at site x. Although we refer to these particles as 
"charged" , they do not interact via the Coulomb potential. Instead, they respond to an external, uniform "electric" 
field E, directed along a specific lattice axis which we label the +xu -direction. We restrict ourselves to zero total 
charge, i.e., X) x [ n x — n x] = 0- I n the absence of the drive, the dynamics does not distinguish between the different 
species: both types hop randomly to nearest neighbor empty sites, with a rate T. The electric field introduces a bias 
into these rates in such a way that jumps against the force will be exponentially suppressed. During one Monte Carlo 
step 11? nearest neighbor bonds are chosen randomly. If a particle-hole pair is encountered, an exchange takes place 
with probability W = T min{l, exp(qE 5x\\ )}, where q = ±1 is the charge of the particle and <5x|| = ±1, is the change 
of the X|| coordinate of the particle due to the jump. 

For our simulations, we set T = 1. Using 60 x 60 lattices, we initialize the system with random configurations of 
various particle densities and carry out runs ranging from 2.5 to 5 x 10 5 MCS. After allowing 62500 MCS for the 
system to settle into a steady state, we measure the Fourier transforms of every 125 MCS, defining them in the 
usual way: 

X 

In the literature, the term "structure factor" typically refers to the (ensemble- or time-) average of density-density 
operators in momentum space, i.e., 

S^(k) = ^K/ k ) , (2) 

where a, /? = +,— ; k = ^(m^my) and V = L 2 is the volume. Instead of the overall averages, we construct 

histograms lor v , ^ — — and ^ — — from their time series m the steady state. We will use the somewhat 

loose term "structure factors" for these fluctuating quantities themselves. By symmetry, k y k and k y ~ k are 
distributed identically, and we consider only the former. The results are presented in Figure [j] and [2] for the smallest 
longitudinal and transverse wave vectors. Before discussing the data, we will first present the theoretical framework 
within which they can be understood. 



III. COARSE-GRAINED DESCRIPTION 



Being interested in the large distance (small k) behavior, we may exploit the simplest approach, i.e., to use coarse- 
grained equations of motions for the slowly varying local densities £» ± (x, t). These can be considered as the naive 
continuum limits of the average lattice occupation numbers (nj) [RH]- For generality, we consider the d-dimensional 
case when X|| is directed along the electric field and xx lies in the d — 1 dimensional subspace transverse to the field. 
Then the mean-field equations of motion for the charge densities follow as: 

dtg* = -VT{g ± V (1 - q+ - q-) ± £X||g ± (l -g+ - g-)}; (3) 

Here, T is the diffusion matrix which is diagonal but not a simple multiple of the unit matrix, due to the anisotropics 

induced by the bias. V is the asymmetric gradient operator, acting on any two functions / and g according to 

/ V g = /Vg — <?V/. e is the coarse-grained bias and xy is the unit vector along the x\\ direction. The coefficients 
appearing in Eqn. (|3|) are functions of the microscopic control parameters; however, their precise dependence need 
not be known. 

To study small fluctuations in the disordered phase, we linearize (||) about the homogeneous time-independent 
solution g and add conserved noise terms, to model the fast degrees of freedom. Writing g ± (x, t) = g + x ± ( x , t), we 
arrive at the Langevin equations for the density fluctuations: 

9 tX a (x,t) = /: a/3 (V)x' 3 (x,t)-Vr 7 Q (x,i) , (4) 

where 
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rz^vii- Ai-e)vrv-(i-3e)er||0u evrv + 0er|^| \ 

1 1 jj_ V gvrv-per||9|| (i - e)vrv + (l - 3^r||an j lj 

and summation over repeated indices is implied in (^|) and in the following. Here, ?7 ± (x, t) are Gaussian white noise 
terms, satisfying: 

(ijf(x,i)) =0, 

«(x, t)^(x', f )> = 5(x - x')*(t - f ) , W 

where a, /3 = +, — ; i, j = 1,2, . . . d and £ Q/3 is the Kronecker-symbol. (<Xy t) = <x is the noise matrix, and similar to T, 
it is diagonal but not proportional to the unit matrix, due to the bias. To describe systems in thermal equilibrium, 
the fluctuation-dissipation theorem would require <x oc T. In particular, in the absence of the drive, we would have 
er = p{l — 2p)r here. However, when driven, the proportionality is not expected to hold, since the diffusion and noise 
matrices are renormalized differently by the drive e, similar to the driven single species case Q. Finally, we point out 
that r} + and r\~ are uncorrelated, due to the fact that charge exchange is not allowed. 

Introducing the Fourier transforms for the density fluctuations (and similarly for the noise) 

X ± (k, u) = J dtd d x X ± (x, f ) e -*M+ kx ) , (7) 

the solution to (0) is trivial: 

X Q (k,c) = (i- 1 )^ l V(k, W ), (8) 

where L a @ — C af3 (ik) — iu5 a ^ . According to the definition of the full dynamic structure factors 

5^(k, u) [(2vr) d+1 ( 5(k + V)6(w + w ')] = ( x «(k, cj) X ?(k', w ')> , (9) 

we can integrate them over u> to obtain the steady state structure factors: 

S^(k) [(27r) d ,5(k + k')] = (x"(k, t) X P {k', t)) . (10) 

Note that they are just the Fourier transforms of the usual equal-time correlation functions (x"(x, t)x (0> t)) = 
G"^(x) . The explicit forms of the structure factors are: 

, ++(k) = izi k^k (krk) 2 + i^f rf 
[ ' (i - 2§) krk (krk) 2 + (i - 4e)rjj£ 2 fc2 

kak (krk) 2 - |g| r 2 £ 2 fcf 

(1 - 2q) krk (krk) 2 + (I - 4g)r 2 e 2 /c 2 



Re{S+-(k) } = - < = ^ - " :» (11) 



Im{5 + -(k)} 



)2 _)_ n _ /l^p2 ff 2/,2 

2g k<Tkri|£fc|i 



(i - q) (krk) 2 + (i - 4g)r 2 £ 2 fc 2 ' 



where k= (kj_,fe||). They clearly exhibit a finite discontinuity singularity at the origin, i.e., limfe,,—^ 5 ,Q/3 (0, k\\) ^ 
limk^^o ^"^(kj,, 0) which translates to power law decays in the spatial correlations. In particular, we find the 
expected fl7j| power law, r~ d , related to er <jt T, typical of non-equilibrium steady states of a system with anisotropy 
and subjected to a conservation law. In addition, a novel power, r~( d+1 ^ 2 controls correlations along the bias 
direction. 

Next, we seek the full probability distributions of the density-density operators. Thus, for each k vector, we construct 
the distributions separately for x ( k,t )* (~ k ' t ) ; Re[x (k,t)x (-M)] an( j Im[ x (k,t)x (-M)] ^ re p resen ted by s ++ , s+~ 
and : 

++. w _ / x / X + (k,<)x + *(k, f) _ +H 



P ++ ( S ++ ; k) = ( S [ * y ' 'f v ' 7 - a 

pr «-;k) = ( * (' Re[x+(k '^^ (k '^ )] - .r i ) (12) 



i?~(*+-;k) = ((S 



Im[x + (k,t)x-*(k,t)] 
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The normalization by V helps to avoid infinities in the expectation values of s ++ , s+ , s+ as can be seen from ( JTc| ) 
and by noting that (2ir) d S(k = 0) — V in the infinite volume limit. Also note that we used x ± (— k, t) = (k, t), since 
the densities x ± (r,t) are real. For our purposes, we need both, the density fluctuations and the explicit distribution 
of the noise, in the (k, t) domain. Within the linear stability regime of the disordered phase, each matrix element of 
{(L~ 1 ) a ^) has two poles in the positive to half-plane corresponding to two stable eigenvalues of (£ a/3 ). Thus, from 
(§) we find: 



x°(M) 



dt' 



2ir 



tkafQt, t') 



(13) 



with 



P[rjf(k,t)] oc exp- 



(14) 



where Df? = 25 a P<Tij in our model. Due to translational invariance, fields with different k vectors are decoupled, so 
we will suppress k in the following. Then (13) can be written as 



X a {t)=vf (t,t')r,?(t>) 



where 



vf(t,t') = Q(t-t') 



(15) 



(16) 



Note that summation over repeated indices also includes an integral over t' in ( |T5| ) 

x + (t)x + '(t) 

V 



We start with the probability distribution of 



by first finding its characteristic function, i.e., 



ds ++ e ins++ - I 



P++(s ++ ) = (e 



(17) 



When performing the average in ([l7j), all integrations over the noise are trivial, except those associated with ±k. 
Thus, we need to evaluate the following integral: 



P ++ (n) = J J] d V ](t')dr,r(t') P[^«(t), ,?*(*)] e fl 



(18) 



where T]^(t') and rjj*(t r ) are the k and — k components of the noise, respectively, yielding the only non-trivial inte- 



grations. Inserting ( |l5| ) into (|lj), we still have a Gaussian integrand, controlled by the quadratic form: 



((D- 1 ) flv -inv+v+*)i lv * 



(19) 



where the indices /i, v include all the degrees of freedom left over, namely, time, charge, and spatial component. The 
corresponding path integrals lead to 



P++(Q) = 



det [(D- 1 )^] 



1 



det 



det 



5^ - iQD^VyVu 



(20) 



where the numerator in the middle expression originates in the normalization factor ensuring P ++ (0) = 1. Exploiting 
the formula det [8^ v + a^b w ] = 1 + a^b^, we obtain 



1 



1 - iflvpDnvVv 



Note that the coefficient of zO is simply the '++' structure factor: 

'x + (M)x + *(k,t) 



V 



s++(k) 



(21) 



(22) 
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Taking the inverse transform to obtain P ++ , the single pole — i/S ++ (k) in the lower half fl— plane yields an exponential 
distribution for the non- negative variable s ++ , i.e., 



4nrr e~ s++ / S++ ^ if S++ > 
if s++ < 



P++(.s++;k) = { sT+(k) e ^ u *^ ^ u (23) 



We will refer to l/S' ++ (k) as the "decay factor" of the exponential. 

Next, we consider the distribution of "^ e ^ x y£ — — . Following the same steps, we see that the quadratic form 
controlling the Gaussian integrand is 



(24) 



leading to 



P<- (n) = — r ~ ; — : ~rr • ( 25 ) 



det 



Using the formula det [S^ v + a^b v + c^d,,] = 1 + [a^b^ + c^d^} + [{a^b^){c v d u ) - (a Al d M )(6 l/ c l/ )] we arrive at 

1 

l-/OR>.,s ik)] !-if[|S (kil-'-l.s^ (k) -'] 
where, similarly to (f22), we have used 



Pr "(»)= , .»„ ,^_„„ , a ,„„^„ ,„ ,^_„»„, ■ (26) 



,; v ,;.(>ifM).,- w . (27) 

Unlike the previous case, P r + (0) has two poles: one (0_) being on the negative, and the other (f2+) on the positive, 
imaginary axis: 

n=F = | (Re[5+"(k)] =F VA+(fe[S+-(k)])2) (28) 

where A = |5 ++ (k)| 2 — (k)| 2 > . In general, their magnitudes are different. Thus, the inverse transform yields 
an asymmetric exponential distribution, characterized by two distinct decay factors and 



P+-(,+-;k) = AT e ; S r_ £ U (2Q) 

I .A. e l"+|S r if SJ < 



iV 



with TV = + (Rc[5+-(k)]) 2 . To obtain the distribution of Im[x+( ^ we simply interchange Rc[5 +_ (k)] and 
Im[5+-(k)] in cqs. (||||). 

One important consequence of ( p3f ) and (29) is that their standard deviations always take a value greater or equal 
than the average. In particular 



^/(a ++ ) 2 -(*++) 



2 = S++(k) 



(>" )" = \/(Re[5+-(k)]) 2 + ^ (so) 



^-) 2 -(^-) 2 = ^/(Im[5+-(k)])2 



A 



Thus, when measuring structure factors in the disordered phase, fluctuations comparable to the average should not 
come as a surprise. We also emphasize that formulas (|||) and (|2^) are completely independent of the specific model 
and can be derived in any multicomponent system. The only necessary conditions are linear Langevin equations with 
Gaussian noise. 
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IV. DISCUSSION AND SUMMARY 



Finally, let us turn to a comparison with the simulation results summarized in Figure [l] and [2| for the two smallest 
wave vectors. The '++' histograms show the usual exponential decay ||], while the 'H — ' histograms clearly represent 
asymmetric exponential distributions. To test the predictions of our Gaussian theory, namely, that the slopes of the 
histograms are determined by the structure factor averages themselves, we simply measured these averages, i.e., S ++ 
, ReS^ and ImS* -1 . The '++' case is particularly simple since the decay factor is just the inverse of S ++ itself. For 
the two 'H — ' distributions, we inserted the measured averages into the non-trivial theoretical relationship (p8|), to find 
the decay factors |f2=p|. The agreement between simulation results and theory is remarkably good, indicating that 
linearized Langevin equations are quite acceptable in this regime. While the external field, E, may obviously generate 
renormalizations, these can be absorbed in the effective parameters of the theory, namely the diffusion matrix, the 
noise matrix, the average density and the coarse-grained bias, leaving the form of the structure factor distributions 
invariant. However, we must avoid critical fluctuations since the linear approximation will fail to capture their effects 
correctly. 

To summarize, we have analyzed the full structure factor distributions in a simple driven lattice gas. Due to the 
presence of two species, there are several distinct structure factor distributions here. All are simple exponentials, with 
the decay factors given in terms of the average structure factors. However, only the "diagonal" elements P ++ and 
P are characterized by a single decay factor, being distributions for a non-negative variable. In contrast, both Pp 
and P+- are combinations of two exponentials, with different decay factors, characterizing positive versus negative 
arguments, respectively. While this asymmetry may, at first sight, be rather surprising, its origin is quite easily traced 
to the presence of off-diagonal elements in the matrix £ a @, Eqn. (|^). Thus, the observation of asymmetric structure 
factor histograms, either in an experiment or a simulation, provides some basic insight into the form of the underlying 
Langevin equation. 
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FIG. 1. Histograms representing the distributions of the k = ^(0, 1) structure factors for (a) k v k , (b) ^ — — and 

(c) Im[ " k "- kl . L = 60, E = 0.471 and g = 0.175. Theoretical distributions (a) P++(s ++ ;k), (b) P+"(s+-;k) and (c) 
P ; H (s^ ; k) are plotted with solid lines on the same graphs. 
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FIG. 2. Histograms representing the distributions of the k = 2e(1,0) structure factors for (a) k v k , (b) ^ — — and 

(c) Im[ " k "- kl . L = 60, E = 0.471 and g = 0.175. Theoretical distributions (a) P++(s ++ ;k), (b) P+"(s+-;k) and (c) 
P ; H (s^ ; k) are plotted with solid lines on the same graphs. 
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